Method for the Reconstruction of Compressed Digital Image Data

ABSTRACT

The invention relates to a method for the re-construction of digital image data which has been lossy-compressed by transform coding and quantization, comprising the steps of a) Inputting JPEG-compressed digital image data (J comp ); b) Iteratively performing updates of an extended image model, thereby reducing the total generalized variation functional constrained to all images matching the transform coded, quantized data; c) Outputting a decompressed digital image (u).

FIELD OF THE INVENTION AND DESCRIPTION OF PRIOR ART

The invention relates to a method for the reconstruction of compressed digital image data which has been lossy-compressed by transform coding and quantization.

With the expanding use of data communication and increasing computer processing power the need for efficient methods to store and transmit image data is omnipresent. In order to do so image compression is of paramount importance.

A fundamental goal of image compression is to reduce the amount of data necessary to represent an image while maintaining acceptable image quality and correctly reconstructing the original (i.e., uncompressed) image. On the other hand, it is important to compress and decompress image data using few international standards so that communication and transfer of image data between different operating systems and equipment solutions is possible.

There are lots of different standards available today, e.g. JPEG (Joint Photographic Experts Group), JPEG2000, JBIG (Joint Bilevel Image Group), GIF (Graphics Interchange Format) and PNG (Portable Network Graphics), to name only a few. The following considerations are based on the JPEG-standard because it is the most commonly used method for compression in the digital domain. However, it has to be noted that the solution according to the invention as described below is applicable to different compression and decompression methods.

As an introduction to the area of image decompression the JPEG-standard ISO/IEC 10918-1:1994 is explained briefly as an example: FIG. 1A shows an exemplary flow-chart of how a digital uncompressed representation of an original image 100 is transformed into a representation 101 that allows identifying and excluding data that is not necessary for visual perception like high frequency variations of colour and intensity. The numbers in FIG. 1A are arbitrary and serve to give a rough idea of what happens to image data while JPEG-compression.

The transformation is achieved by quantising (step 102) the data and subsequently rounding (step 103) it to integers. Since many possible real numbers are rounded to the same integer the rounded compressed data 104 does not represent the one original image 100 but a set of possible original images. The resulting information is then compressed and stored without further losses.

FIG. 1B shows a supplemental block diagram of abovementioned JPEG-compression. Source image data 100 is compressed by firstly applying a DCT (Discrete Cosine Transformation) on 8×8-blocks of the data. Subsequently, the data is quantized 102, using a quantization table 102′, and rounded 103, followed by a lossless compression 103′ resulting in the lossy-compressed transform coded quantized data 104.

By direct inversion of this process standard reconstruction techniques reconstruct only the image that would directly lead to integer coefficients anyway. Said techniques do not take into account that the inverted data results from a rounding process.

Hence, the JPEG-standard for digital images (ISO/ IEC 10918-1:1994) is lossy, that is, not the entire original image is stored but only the data that is critical to the visual impression. This allows for a higher compression rate but leads to errors—so called artifacts—in the reconstructed image negatively affecting the image quality when conventional reconstruction techniques are used.

In order to compensate for abovementioned errors conventional reconstruction techniques apply postprocessing measures—by that means quality may be increased but the reconstructed image differs from the original image because postprocessing may change the data.

A solution that accounts for this problem is presented in U.S. Pat. No. 5,534,925: Here, it is acknowledged that a JPEG-compressed image comprises many images that are consistent with the given information. Hence, the reconstructed image is selected from all images matching the data as the image that minimizes the sum of the difference of each pixel to its neighbouring pixels. In the context of a transform-quantization coding it is stated that the optimal reconstruction is a nonlinear optimization with linear inequality constraints. The most successful criterion used for image processing is minimizing oscillation—minimizing the total variation (TV) with constraints minimizes oscillations while maintaining the edges of structures in the image sharp. This, however, leads to stairstepping-effects (aliasing-effects) negatively influencing image quality. Hence, it is suggested to cancel the algorithm after a few iterations.

It is well accepted that methods based on TV minimization lead to abovementioned stairstepping- or staircasing effexts, resulting in unnatural reconstructed images. If the image to be reconstructed contains not only flat but also slanted regions, then the image reconstructed on the basis of the bounded variation seminorm tends to be piecewise constant (staircasing). No methods solely based on or similar to TV-minimization and the resulting “classical” image model, where only brightness (or color) information of an image are considered, are able to overcome abovementioned drawbacks without introducing new disadvantages.

SUMMARY OF THE INVENTION

The invention sets out to overcome the above-mentioned shortcomings of the prior art by providing a method that reconstructs lossy-decompressed image data so that it is close to the original data with a minimum amount of errors.

This task is solved by a method according to the invention, comprising the steps of

-   -   a) Inputting a lossy-compressed digital image data;     -   b) Iteratively performing updates of an extended image model,         thereby reducing the total generalized variation functional         constrained to all images matching the transform coded,         quantized data;     -   c) Outputting a decompressed digital image.

By virtue of this solution it is possible not to amend a reconstructed image in order to enhance its quality, but merely to choose the optimal image from the set of possible images—the method accounts for the fact that a JPEG-compressed image may emanate from a large set of original images and that this information is lost if not considered during reconstruction.

The total generalized variation model is discussed, e.g., in the article “Total Generalized Variation” by Bredies et al., Society for Industrial and applied Mathematics—J. Imaging Sciences, 2010, Vol. 3, No. 3, pp. 492-526. The disclosure of said article is enclosed to this specification by reference.

The invention relates to a method for optimal and artifact-free reconstruction of digital image data compressed by transform coding and quantization. The image data may come from a digital camera, e.g. in a mobile phone, but also application in web browsers or video streams to increase image quality is possible. The reproduction process iteratively performs updates of an extended image model, represented by primal and dual variables, respectively, thereby reducing the total generalized variation functional constrained to all images matching the transform coded, quantized data. Coefficient and quantization data allow description of a coefficient data set and, consequently, an image data set.

Using the TGV functional means including not only brightness (or color), but also gradient information of arbitrary order as part of the image and allows an optimal balancing between the first- and higher order derivatives. This is the fundamental difference to other methods, such as TV-based methods, which only measure first order derivatives. As a consequence, the method according to the invention significantly improves image quality compared to known reconstructon models, e.g., by reducing the staircasing effect of the bounded variation functional.

The usage of an extended image model results from the nature of the TGV functional. Depending on the order of the TGV functional in its most general form it comprises brightness (or color) and gradient information up to a certain order. A variant to perform TGV minimization using this extended image model is to include abovementioned brightness and gradient information in primal and dual variables. Aoso other variants are possible, however, usage of the extended image model in some form is implied by TGV minimization and thus essential fore the highly improved reconstruction quality.

The method according to the invention does not post process decompressed data but reconstructs an image fitting the compressed data. By application of the novel TGV-model an image is selected automatically that comes very close to the original image.

The method allows for a quick and optimal reconstruction of a compressed JPEG-image with the reconstruction being guaranteed by the mathematical validity of the TGV-model and by convex analysis. The method does not have the drawbacks of existing methods. It results in natural and visually more appealing image reconstructions.

In step a) coefficient data (d) and quantization data (Q) is obtained from the lossy transform-coded digital image data (J_(comp)) and the coefficient data (d) is dequantized. This is, for instance, a standard procedure when opening a JPEG-compressed image. The quantization data (Q) usually is provided in the form of a table.

Coefficient data here denotes the quantized DCT (Discrete cosine transformation)-coefficient matrices for each color component or the brightness component in case of color- and grey-scale images, respectively. In case of a color image the matrices can be seen as three scalar-valued matrices corresponding to greyscale data. The quantization data is a composition of quantization matrices corresponding to the three color components or the brightness component.

In a variant of the invention, in step b) the extended image model comprises an image, an extended image, a first dual variable and a second dual variable, an extragradient image and an extragradient extended image.

Coefficient and quantization data allow description of a coefficient data set and, consequently, an image data set.

In a variant of the invention, step b) comprises the following sub-steps:

-   -   b0) Initializing image data (u), extended image data (v),         extragradient image data (ū), extragradient extended image data         ( v), a two-vector valued matrix (p) and a three-vector valued         matrix (q) which are used as a starting point for the iteration         procedure in steps b1) to b8);     -   b1) Applying a gradient evaluation and a         summation/multiplication operation on the extragradient image         data (ū), the extragradient extended image data ( v) and the         two-vector valued matrix (p), thereby updating the two-vector         valued matrix (p) and then applying a projection operation on         said matrix (p) in a way that its norm (n1_(i,j)) is less or         equal to a first parameter (α₁) and each two-vector entry of the         two-vector valued matrix (p) is modified accordingly;     -   b2) Applying a higher dimensional gradient evaluation and a         summation/multiplication operation on the extragradient extended         image data ( v) and the three-vector valued matrix (q), thereby         updating the results in the three-vector valued matrix (q) and         then applying a projection operation on said matrix (q) in a way         that its norm (n2_(i,j)) is less or equal to a second parameter         (α₀) and each three-vector entry of the three-vector valued         matrix (q) is modified accordingly;     -   b3) Applying a divergence operator and a         summation/multiplication operation on the image data (u) and the         two-vector valued matrix (p) of step b1), storing the data in a         first new set of image data (u_(new));     -   b4) Applying a divergence operator and a         summation/multiplication operation on the extended image data         (v), the two-vector valued matrix (p) of step b1) and the         three-vector-valued matrix (q) of step b2), storing the data in         a second new set of image data (v_(new));     -   b5) Updating the first new set of image data (u_(new)) pf step         b3) by applying a sub-sampling operation S, a component-wise         block-wise discrete cosine transformation and a projection         operation on the first new set of image data, followed by an         inverse component-wise block-cosine transformation IBDCT, an         up-sampling operation S⁻¹ and a summation operation;     -   b6) Updating the extragradient image data (ū) and the         extragradient extended image data ( v) by an         extragradient-update step using the image data (u), extended         image data (v), first (u_(new)) and second new set of image data         (v_(new));     -   b7) Overwriting image data (u) and extended image data (v) with         of the first (u_(new)) and second new set of image data         (v_(new)), respectively, resulting from the steps b1) to b6);     -   b8) Evaluating a stopping-criterion;     -   b9) Reiterating steps b1) to b8) if stopping criterion is not         fulfilled or proceeding to step c) if stopping criterion is         fulfilled.

The two-vector valued matrix (p) and the three-vector valued matrix (q) are abovementioned first and second dual variable, respectively.

In step b0) different types of starting images may be used. In principle, any initialization is suitable, as for example also an image {right arrow over (u)}=0, however, in a variant of the invention a standard-decompressed image is selected, performing an inverse component-wise block-cosine transformation IBDCT and an up-sampling operation S⁻¹ on the dequantized data (d).

In a variant of the invention in step b1) the gradient-evaluation maps extragradient image data ( {right arrow over (u)}) to a two-vector-valued matrix with the entries being defined as follows:

( V {right arrow over (u)}) _(i,j) ^(n,1)=(δ_(x+) ū ^(n))_(i,j), and

( V {right arrow over (u)}) _(i,j) ^(n,2)=(δ_(y+) ū ^(n))_(i,j),

with 0≦i<8 k and 0≦j<8 l, k and l corresponding to the vertical and horizontal number of 8×8 data blocks and n=1, 2, 3, the operators each resulting in a scalar-valued matrix with

$\left( {\delta_{x +}\overset{\_}{u}} \right)_{i,j} = \left\{ {{\begin{matrix} \left( {{\overset{\_}{u}}_{{i + 1},j} - {\overset{\_}{u}}_{i,j}} \right) & {{{{if}\mspace{14mu} 0} \leq i < {{8k} - 1}},} \\ 0 & {{{{if}\mspace{14mu} i} = {{8k} - 1}},} \end{matrix}{{and}\left( {\delta_{y +}\overset{\_}{u}} \right)}_{i,j}} = \left\{ {\begin{matrix} \left( {{\overset{\_}{u}}_{i,{j + 1}} - {\overset{\_}{u}}_{i,j}} \right) & {{{{if}\mspace{14mu} 0} \leq j < {{8l} - 1}},} \\ 0 & {{{if}\mspace{14mu} j} = {{8l} - 1}} \end{matrix}{and}} \right.} \right.$

the summation/multiplication operation yields a two-vector valued matrix (p) with

p _(i,j) ^(m,n) =p _(i,j) ^(m,n)+σ(( V {right arrow over (u)}) _(i,j) ^(m,n) − v _(i,j) ^(m,n)),

σ being a first scalar and v _(i,j) ^(m,n) being the extragradient extended image data ( {right arrow over (v)}).

In yet another variant of the invention, in step b1) the modification of the two-vector valued matrix (p_(i,j) ^(m,n)) is calculated with

${p_{i,j}^{m,n} = \frac{p_{i,j}^{m,n}}{\max \left\{ {1,\frac{n\; 1_{i,j}}{\alpha_{1}}} \right\}}},$

α₁ being a first parameter and n1 being a norm.

The norm (n1_(i,j)) is calculated as

${n\; 1_{i,j}} = {\sqrt{\left( p_{i,j}^{1,1} \right)^{2} + \left( p_{i,j}^{2,1} \right)^{2} + \left( p_{i,j}^{3,1} \right)^{2} + \left( p_{i,j}^{1,2} \right)^{2} + \left( p_{i,j}^{2,2} \right)^{2} + \left( p_{i,j}^{3,2} \right)^{2}}.}$

This is one possible norm given by the TGV-model.

In step b2) the higher dimensional gradient evaluation maps extragradient extended image data ( v) to a three-vector valued matrix being defined as follows:

${\left( {ɛ\; \overset{\overset{\rightarrow}{\_}}{v}} \right)_{i,j} = \left( {\begin{pmatrix} \left( {\delta_{x -}{\overset{\_}{v}}^{1,1}} \right)_{i,j} \\ \left( {\delta_{x -}{\overset{\_}{v}}^{2,1}} \right)_{i,j} \\ \left( {\delta_{x -}{\overset{\_}{v}}^{3,1}} \right)_{i,j} \end{pmatrix},\begin{pmatrix} \left( {\delta_{y -}{\overset{\_}{v}}^{1,2}} \right)_{i,j} \\ \left( {\delta_{y -}{\overset{\_}{v}}^{2,2}} \right)_{i,j} \\ \left( {\delta_{y -}{\overset{\_}{v}}^{3,2}} \right)_{i,j} \end{pmatrix},\begin{pmatrix} \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{1,1}} + {\delta_{x -}{\overset{\_}{v}}^{1,2}}}{2} \right)_{i,j} \\ \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{2,1}} + {\delta_{x -}{\overset{\_}{v}}^{2,2}}}{2} \right)_{i,j} \\ \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{3,1}} + {\delta_{x -}{\overset{\_}{v}}^{3,2}}}{2}\; \right)_{i,j} \end{pmatrix}} \right)},$

the operators being defined as

$\left( {\delta_{x -}\overset{\_}{v}} \right)_{i,j} = \left\{ {{\begin{matrix} {- {\overset{\_}{v}}_{{i - 1},j}} & {{{if}\mspace{14mu} i} = {{8k} - 1}} \\ \left( {{\overset{\_}{v}}_{i,j} - {\overset{\_}{v}}_{{i - 1},j}} \right) & {{{if}\mspace{14mu} 0} < i < {{8k} - 1}} \\ {\overset{\_}{v}}_{i,j} & {{{if}\mspace{14mu} i} = 0} \end{matrix}{{and}\left( {\delta_{y -}\overset{\_}{v}} \right)}_{i,j}} = \left\{ \begin{matrix} {{- {\overset{\_}{v}}_{i,j}} - 1} & {{{if}\mspace{14mu} i} = {{8l} - 1}} \\ \left( {{\overset{\_}{v}}_{i,j} - {\overset{\_}{v}}_{i,{j - 1}}} \right) & {{{{if}\mspace{14mu} 0} < j < {{8l} - 1}},} \\ {\overset{\_}{v}}_{i,j} & {{{if}\mspace{14mu} j} = 0.} \end{matrix} \right.} \right.$

and the summation/multiplication operation yields a three-vector valued matrix (q_(i,j) ^(m,n)) with q_(i,j) ^(m,n)=q_(i,j) ^(m,n)+σ(ε {right arrow over (v)})_(i,j) ^(m,n).

In step b2) the the modification of the three-vector valued matrix (q_(i,j) ^(m,n)) is calculated with

${q_{i,j}^{m,n} = \frac{q_{i,j}^{m,n}}{\max \left\{ {1,\frac{n\; 2_{i,j}}{\alpha_{0}}} \right\}}},$

α_(o) being a second parameter and n2 being a norm. Here, the norm (n2_(i,j)) is calculated with

${n\; 2_{i,j}} = {\sqrt{\begin{matrix} \begin{matrix} {\left( q_{i,j}^{1,1} \right)^{2} + \left( q_{i,j}^{2,1} \right)^{2} + \left( q_{i,j}^{3,1} \right)^{2} +} \\ {\left( q_{i,j}^{1,2} \right)^{2} + \left( q_{i,j}^{2,2} \right)^{2} + \left( q_{i,j}^{3,2} \right)^{2} +} \end{matrix} \\ {2\left( {\left( q_{i,j}^{1,3} \right)^{2} + \left( q_{i,j}^{2,3} \right)^{2} + \left( q_{i,j}^{3,3} \right)^{2}} \right)} \end{matrix}}.}$

This is yet another possible norm given by the TGV-model.

The first and second parameter α₀ and α₁ have to fulfil the condition α₀, α₁>0.

In a variant of the invention in step b3) a divergence operator is applied to the two-vector valued matrix (p_(i,j) ^(m,n)) of step b1) as follows:

(div {right arrow over (p)})_(i,j) ^(n)=(δ_(x−) p ^(n,1))_(i,j)+(δ_(y−) p ^(n,2))_(i,j), and the summation/multiplication operation yields a first new set of image data ({right arrow over (u)}_(new)) as a vector-valued matrix with

(u _(new))_(i,j) ^(n) =u _(i,j) ^(n)+τ(div {right arrow over (p)})_(i,j) ^(n),

τ being a second scalar.

Furthermore, in step b4) a divergence operator is applied to the three-vector valued matrix (q_(i,j) ^(m,n)) of step b2) as follows:

(div {right arrow over (q)})_(i,j) ^(n,1)=(δ_(x−) q _(n1))_(i,j)+(δ_(y+) q ^(n,3))_(i,j),

(div {right arrow over (q)})_(i,j) ^(n,2)=(δ_(x+) q ^(n,3))_(i,j)+(δ_(y+) q ^(n,2))_(i,j)

and the summation/multiplication operation yields a second new set of image data ({right arrow over (v)}_(new)) as a two-vector valued matrix with

(v _(new))_(i,j) ^(m,n) =v _(i,j) ^(m,n)+τ({right arrow over (p)}+div {right arrow over (q)})_(i,j) ^(m,n).

For abovementioned variables it holds that the first (σ) and the second scalar (τ) satisfy the condition στ≦ 1/12.

In a variant of the invention, in step b5) a sub-sampling operation is applied to the first new set of image data ({right arrow over (u)}_(new)) resulting in scalar-valued matrices as follows:

${\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{1} = {\frac{k_{1}l_{1}}{kl}{\sum\limits_{m = {i\; \frac{k}{k_{1}}}}^{{{({i + 1})}\frac{k}{k_{1}}} - 1}{\sum\limits_{n = {j\; \frac{l}{l_{1}}}}^{{{({j + 1})}\; \frac{l}{l_{1}}} - 1}\left( u_{new} \right)_{n,m}^{1}}}}},,{\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{2} = {\frac{k_{2}l_{2}}{kl}{\sum\limits_{m = {i\; \frac{k}{k_{2}}}}^{{{({i + 1})}\; \frac{k}{k_{2}}} - 1}{\sum\limits_{n = {j\; \frac{l}{l_{2}}}}^{{{({j + 1})}\; \frac{l}{l_{2}}} - 1}\left( u_{new} \right)_{n,m}^{2}}}}},{\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{3} = {\frac{k_{3}l_{3}}{kl}{\sum\limits_{m = {i\; \frac{k}{k_{3}}}}^{{{({i + 1})}\; \frac{k\;}{k_{3}}} - 1}{\sum\limits_{n = {j\; \frac{l}{l_{3}}}}^{{{({j + 1})}\; \frac{l}{l_{3}}} - 1}\left( u_{new} \right)_{n,m}^{3}}}}}$

with k₁, h, k₂, k₃, l₂, l₃ being the sub-sampled image size parameters, then a component-wise block-wise discrete cosine transformation is performed resulting in scalar-valued matrices (BDCT(S{right arrow over (u)}_(new)))¹, (BDCT(S{right arrow over (u)}_(new)))² and (BDCT(S{right arrow over (u)}_(new)))³ which are then processed in a projection operation with

${{proj}_{d,Q}\left( {{BDCT}\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)} \right)}_{i,j}^{n} = {{{BDCT}\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)}_{i,j}^{n} - {\max {\left\{ {0,{{{BDCT}\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)}_{i,j}^{n} - {\left( {d_{i,j}^{n} + \frac{1}{2}} \right)Q_{i,j}^{n}}}} \right\}++}\max \left\{ {0,{{\left( {d_{i,j}^{n} - \frac{1}{2}} \right)Q_{i,j}^{n}} - {{BDCT}\left( {S\; {\overset{\rightarrow}{u}}_{new}} \right)}_{i,j}^{n}}} \right\} \mspace{14mu} {for}}}$

n=1, 2, 3, resulting in three scalar-valued matrices which are then inversely transformed IBDCT and inversely sub-sampled S⁻¹, resulting in

proj_(D)({right arrow over (u)}_(new))_(i,j) ^(n)=(u_(new))_(i,j) ^(n)+(S ¹(IBDCT(proj_(d,Q)((S{right arrow over (u)} _(new))))−S{right arrow over (u)} _(new)))_(i,j) ^(n)

which is then stored as first new set of image ({right arrow over (u)}_(new)).

The extragradient update-step is performed by processing image data ({right arrow over (u)}, {right arrow over (v)}) and the first ({right arrow over (u)}_(new)) and second new set of image data ({right arrow over (v)}_(new)) as follows:

ū _(i,j) ^(n)=2(u _(new))_(i,j) ^(n) −u _(i,j) ^(n),

v _(i,j) ^(m,n)=2(v _(new))_(i,j) ^(m,n) −v _(i,j) ^(m,n)

resulting in new sets of extragradient image data ( {right arrow over (u)}) and extragradient extended image data ( {right arrow over (v)}).

The method according to the invention is applicable on different types of images. The JPEG-compressed digital image data (J_(comp)) may be

-   -   a greyscale image consisting of a matrix (u_(i,j)), each matrix         element indicating the intensity of the respective pixel, or     -   a color image in an RGB-color-space-setting consisting of a         three-vector valued matrix

$\left( {\left( {\overset{\rightarrow}{u}}_{i,j} \right) = \begin{pmatrix} u_{i,j}^{1} \\ u_{i,j}^{2} \\ u_{i,j}^{3} \end{pmatrix}} \right),$

-   -   each matrix vector indicating the amount of Red, Green and Blue         light in the corresponding (i, j)-th pixel of the image, or     -   a color image in an YCbCr-color-snace-setting consisting of a         three-vector valued matrix

$\left( {\left( {\overset{\rightarrow}{u}}_{i,j} \right) = \begin{pmatrix} u_{i,j}^{1} \\ u_{i,j}^{2} \\ u_{i,j}^{3} \end{pmatrix}} \right),$

-   -   the matrix vectors indicating the brightness, the         blue-difference and the red-difference of the corresponding         pixel.

BRIEF DESCRIPTION OF THE DRAWINGS

In the following, the present invention is described in more detail with reference to the drawings, which show:

FIG. 1Aa flow-chart of conventional JPEG-compression (prior art) according to ISO/IEC 10918-1:1994;

FIG. 1B a block-diagramm of conventional JPEG-compression (prior art) of FIG. 1A, and

FIG. 2 a flow-chart of the method according to the invention.

DETAILED DESCRIPTION OF THE INVENTION

It should be appreciated that the invention is not restricted to the following embodiment which merely represents one of the possible implementations of the invention. Furthermore, it is noted that the representations in the figures are only schematic for the sake of simplicity.

FIG. 2 shows a flow-chart of the method according to the invention. An arbitrary camera device 105 captures an image which is subsequently JPEG-compressed including the steps of sub-sampling, forward transforming, quantization and encoding. It has to be noted that JPEG-compression is only one of many possible compression methods—the method is not limited to JPEG, this standard is merely chosen here because it is most commonly used. The digital image data is then stored in an adequate memory 106. The reconstruction method according to the invention may be performed in a dedicated reconstruction processor 107 as indicated in FIG. 2 by the dashed line.

A suited algorithm for the implementation of the method according to the invention is based on a primal-dual algorithm as presented by Chambolle and Pock in “A first-order primal-dual algorithm for convex problems with applications in imaging”, Journal of Mathematical Imaging and Vision, Vol 40, Issue 1, pp 120-145.

The method according to the invention allows for optimal and artifact-free reconstruction of digital image data compressed by transform coding and quantization. The reproduction process iteratively performs updates of an extended image model, represented by primal and dual variables, respectively, thereby reducing the total generalized variation functional constrained to all images matching the transform coded, quantized data.

The process uses the primal variable to store pixel wise image information as well as an approximation of the difference of neighbouring pixels, the proximal difference image. It uses the dual variable to store dual information for the proximal difference image as well as dual information for the differences of the neighbouring image, the higher-order proximal difference image. Furthermore, pixel wise algebraic operations, difference operations for neighbouring pixels, the coding transform as well as projections are utilized to perform updates of the primal and dual variable.

The method treats the compressed JPEG-image as a set of possible original images and selects a suitable image from said set. It comprises the automated implementation of such a selection process by application of the TGV-model and execution of a suited minimization-algorithm.

Summarizing, the method according to the invention allows reading of a saved JPEG-image, processing of the data of and optimal reconstruction of the data.

This may happen automatically on a personal computer, with the user opening an image file and thereby initializing abovementioned method.

In the following, the method according to the invention is explained in an embodiment where it is applied to a digital color image. The invention, however, is also applicable to a greyscale image. By and large, the difference between the approaches lies in vector-values as opposed to scalar values forming the image data.

A digital color image of the size 8 k×8l pixel is processed as a vector-valued matrix denoted by ({right arrow over (u)}_(i,j))_(0≦i<8 k, 0≦j<8 l). In case of a greyscale image there would be scalar values instead of vectors. The horizontal and vertical pixel numbers are assumed to be multiples of eight times the subsampling factors for color components.

The values k and l are multiples of k₁, k₂, k₃ and l₁, l₂, l₃, respectively, which depend on the subsampling (see below) performed during the compression process of a JPEG-image. Typical values are k₁=k, k₂=k₃=k/2 and l=1 and l₂=l₃=l/2.

This choice of k₁, k₂, k₃ and l₁, l₂, l₃ illustrates the case where the first component is not sub-sampled, but the other components are. Other variants are possible as well with different values for the respective items.

The vector-valued matrix can be visualized by

${\overset{\rightarrow}{u}}_{i,j} = \begin{pmatrix} {\overset{\rightarrow}{u}}_{0,0} & {\overset{\rightarrow}{u}}_{0,1} & {\overset{\rightarrow}{u}}_{0,2} & \ldots & {\overset{\rightarrow}{u}}_{0,{{8l} - 1}} \\ {\overset{\rightarrow}{u}}_{1,0} & {\overset{\rightarrow}{u}}_{1,1} & \; & \; & {\overset{\rightarrow}{u}}_{1,{{8l} - 1}} \\ {\overset{\rightarrow}{u}}_{2,0} & \; & \; & \; & \vdots \\ \vdots & \; & \; & \; & \; \\ {\overset{\rightarrow}{u}}_{{{8k} - 1},0} & \ldots & \; & \; & {{\overset{\rightarrow}{u}}_{{{8k} - 1},{{8l} - 1}}\;} \end{pmatrix}$

where the value of each entry defines the color of the corresponding pixel. In case of a greyscale image each entry of the matrix would be a single scalar value reproducing the intensity of the respective pixel.

Here, each matrix entry {right arrow over (u)}_(i,j), for 0≦i<8 k, 0≦j<8 l, consists of a three-dimensional vector, i.e.

${\overset{->}{u}}_{i,j} = {\begin{pmatrix} u_{i,j}^{1} \\ u_{i,j}^{2} \\ u_{i,j}^{3} \end{pmatrix}.}$

U can be defined as the set of all color images, i.e.

U={({right arrow over (u)} _(i,j))_(0≦i<8 k, 0≦j<8 l)=(u _(i,j) ¹ , u _(i,j) ² , u _(i,j) ³)_(0≦i<8 k, 0≦j<8 l)}.

In a standard RGB color space setting, the three entries u_(i,j) ¹, u_(i,j) ², u_(i,j) ³ represent the amount of Red, Green, and Blue light, respectively, in the corresponding (i, j)-th pixel of the image. However, JPEG-images are typically processed in a YCbCr color space setting.

This means the first line denotes the brightness of the corresponding pixel while the entries of the second and third line describe the blue-difference and the red-difference, respectively.

This is an equivalent color space setting and images given in either RGB- or YCbCr-color space setting can easily be transformed to the other one. Naturally, other color space settings are possible.

The JPEG-standard commonly uses the YCbCr-color space setting because the human eye is less sensitive to color variations than brightness variations. Hence, the color components of a JPEG compressed image (u_(i,j) ², u_(i,j) ³) are typically saved with a lower resolution than the brightness component u_(i,j) ¹, allowing for smaller files. This is called color subsampling and must be taken into account in the decompression process. However, the JPEG standard also allows processing images in the RGB-color space, but then typically no color subsampling is performed because no color channel should be preferred to the others.

As a first step in the reconstruction process data formatting is performed in the reconstruction processor 107. This means that the data necessary for the reconstruction is obtained from the JPEG-compressed image and variables used in the process are initialized.

Besides image data ({right arrow over (u)}_(i,j))_(0≦i<8 k, 0≦j<8 l) two other types of higher dimensional objects are used in the decompression process.

The first type is a two-vector valued matrix denoted by ({right arrow over (v)}_(i,j))_(0≦i<8 k, 0≦j<8 l), where

${\overset{->}{v}}_{i,j} = {\left( {{\overset{->}{v}}_{i,j}^{1},{\overset{->}{v}}_{i,j}^{2},{\overset{->}{v}}_{i,j}^{3}} \right) = {\left( {\begin{pmatrix} v_{i,j}^{1,1} \\ v_{i,j}^{2,1} \\ v_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} v_{i,j}^{1,2} \\ v_{i,j}^{2,2} \\ v_{i,j}^{3,2} \end{pmatrix}} \right).}}$

The second type is a three-vector valued matrix denoted by ({right arrow over (w)}_(i,j))_(0≦i<8 k, 0≦j<8 l), where

${\overset{->}{w}}_{i,j} = {\left( {{\overset{->}{w}}_{i,j}^{1},{\overset{->}{w}}_{i,j}^{2},{\overset{->}{w}}_{i,j}^{3}} \right) = {\left( {\begin{pmatrix} w_{i,j}^{1,1} \\ w_{i,j}^{2,1} \\ w_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} w_{i,j}^{1,2} \\ w_{i,j}^{2,2} \\ w_{i,j}^{3,2} \end{pmatrix},\begin{pmatrix} w_{i,j}^{1,3} \\ w_{i,j}^{2,3} \\ w_{i,j}^{3,3} \end{pmatrix}} \right).}}$

These objects play a role in the internal decompression process and are not crucial for the method according to the invention; however, they are needed to store higher dimensional data resulting from operations on the image.

From a JPEG-compressed object further information can be obtained, namely the following data:

$\overset{->}{d} = {{\begin{pmatrix} d^{1} \\ d^{2} \\ d^{3} \end{pmatrix}\mspace{14mu} {and}\mspace{14mu} \overset{->}{W}} = {\begin{pmatrix} W^{1} \\ W^{2} \\ W^{3} \end{pmatrix}.}}$

d¹, d² and d³ are the quantized DCT (Discrete Cosine Transformation)-coefficient matrices for each color component with d¹=(d_(i,j) ¹)_(0<i<8 k) ₁ _(, 0<j<8 l) ₁ , d²=(d_(i,j) ²)_(0<i<8 k) ₂ _(, 0<j<8 l) ₂ and d³=(d_(i,j) ³)_(0≦i<8 k) ₃ _(, 0≦i<8 l) ₃ . This is equivalent to three grey-scale matrices which correspond to scalar-valued matrices.

As mentioned above the values k and l are multiples of k₁, k₂, k₃ and l₁, l₂, l₃, respectively, which depend on the color subsampling performed during the compression process of a JPEG-image. Typical values are k₁=k, k₂=k₃=k/2 and l₁=l, l₂=l₃=l/2.

The quantization matrix {right arrow over (W)} is a composition of the quantization matrices W¹=(W_(i,j) ¹)_(0≦i−8, 0≦j<8), W²=(W_(i,j) ²)_(0≦i<8, 0≦j<8) and W³=(W_(i,j) ³)_(00≦i<8, 0≦j<8) corresponding to the three color components. The quantization matrices W¹, W² and W³ are uniform for the entire image, i.e., in the compression process, the data corresponding to each color component is split into blocks of size 8×8 and then each of these blocks is divided by the same color-dependent quantization matrix W^(n), 1≦n≦3.

In order to keep the notation as simple as possible these quantization matrices of size 8×8 are extended to the size of the data for the corresponding color component in the following. To illustrate this measure the quantization matrices are now denoted by Q¹=(Q_(i,j) ¹)_(0≦i<8 k) ₁ _(, 0≦j<8 l) ₁ , Q²=(Q_(i,j) ²)_(0≦i<8 k) ₂ _(, 0≦j<8 l) ₂ and Q³=(Q_(i,j) ³)_(0≦i<8 k) ₃ _(, 0≦j<8 l) ₃ where, for example,

$Q^{1} = \begin{pmatrix} W^{1} & W^{1} & \ldots & W^{1} \\ W^{1} & \; & \; & \vdots \\ \vdots & \; & \; & \; \\ W^{1} & \ldots & \; & W^{1} \end{pmatrix}$

and analogously for Q² and Q³.

{right arrow over (Q)} is then defined to be a matrix obtained as a composition of Q¹, Q² and Q³, i.e.,

$\overset{->}{Q} = {\begin{pmatrix} Q^{1} \\ Q^{2} \\ Q^{3} \end{pmatrix}.}$

Using the given data d and Q a coefficient data set D and an image data set U_(D) can be described.

While the sets D and U_(D) are essential for the motivation of the decompression process and also for some steps performed during the process they will not be needed or saved explicitly in the algorithm of the method according to the invention.

Firstly, data intervals J¹=(J_(i,j) ¹)_(0≦i<8 k) ₁ , 0≦j<8 l ₁ , J²=(J_(i,j) ²)_(0≦i<8 k) ₂ _(, 0≦j<8 l) ₂ and J³=(J_(i,j) ³)_(0≦i<8 k) ₃ _(, 0≦j<8 k) ₃ are calculated which are defined by

$\mspace{20mu} {{J_{i,j}^{1} = {{\left\lbrack {{{d_{i,j}^{1}Q_{i,j}^{1}} - \frac{Q_{i,j}^{1}}{2}},{{d_{i,j}^{1}Q_{i,j}^{1}} + \frac{Q_{i,j}^{1}}{2}}} \right\rbrack \mspace{14mu} {for}\mspace{14mu} 0} \leq i < {8\; k_{1}}}},{0 \leq j < {8\; l_{1}}},{J_{i,j}^{2} = {{\left\lbrack {{{d_{i,j}^{2}Q_{i,j}^{2}} - \frac{Q_{i,j}^{2}}{2}},{{d_{i,j}^{2}Q_{i,j}^{2}} + \frac{Q_{i,j}^{2}}{2}}} \right\rbrack \mspace{14mu} {for}\mspace{14mu} 0} \leq i < {8\; k_{2}}}},{0 \leq j < {8\; l_{2}}},{and}}$ $\mspace{20mu} {{J_{i,j}^{3} = {{\left\lbrack {{{d_{i,j}^{3}Q_{i,j}^{3}} - \frac{Q_{i,j}^{3}}{2}},{{d_{i,j}^{3}Q_{i,j}^{3}} + \frac{Q_{i,j}^{3}}{2}}} \right\rbrack \mspace{14mu} {for}\mspace{14mu} 0} \leq i < {8\; k_{3}}}},{0 \leq j < {8\; {l_{3}.}}}}$

The interval J_(i,j) ^(n) is just the interval of all data values z_(i,j) ^(n) that would result in the same data value d_(i,j) ^(n) when divided by the corresponding quantization value Q_(ij) ^(n) and rounded to an integer value.

With that, the source data set D can be defined, i.e. the set of all transformed, sub-sampled images that would result in the coefficient data {right arrow over (d)} when divided by the quantization matrices Q and rounded to an integer value, as

$D = {\left\{ {{\overset{->}{z} = \begin{pmatrix} z^{1} \\ z^{2} \\ z^{3} \end{pmatrix}},{{{such}\mspace{14mu} {that}\mspace{14mu} z_{i,j}^{1}} \in J_{ij}^{1}},{z_{i,j}^{2} \in J_{ij}^{2}},{z_{i,j}^{3} \in {J_{ij}^{3}\mspace{14mu} {for}\mspace{14mu} {all}\mspace{14mu} i}},j} \right\}.}$

This is just the set of all coefficient data {right arrow over (z)} that would produce the same data matrices {right arrow over (d)}, when divided by the quantization matrices {right arrow over (Q)} and rounded to integer.

For the sake of completeness, the complete set of possible source images, i.e., the set of images that would lead to the same data {right arrow over (d)} when JPEG-compression is applied, can now be defined as

U _(D) ={{right arrow over (u)}∈U such that BDCT (S{right arrow over (u)})∈D},

where BDCT(S{right arrow over (u)}) denotes the resulting data matrices when the component-wise, block-wise discrete cosine transformation operator BDCT is applied to the sub-sampled image S{right arrow over (u)} and S{right arrow over (u)} denotes the resulting image after color sub-sampling of the image data {right arrow over (u)}. However, this set U is not directly needed in the TGV-based JPEG-decompression procedure according to the invention.

After abovementioned explanatory remarks about the data {right arrow over (d)} and {right arrow over (Q)} that can be obtained from a given JPEG object the TGV-based JPEG-decompression process according to the invention is described for color images.

The input of the process is a given JPEG compressed object J_(comp), e.g. digital image data from memory 106. The output is a decompressed image u that may be shown on a display device, e.g. a monitor 108.

In a first step, the first “data formatting” box in FIG. 2, coefficient data d and the quantization matrix Q are obtained from the JPEG object and saved. As already mentioned, the data can be seen as composition of three image-type objects, i.e.

${\overset{->}{d} = \begin{pmatrix} d^{1} \\ d^{2} \\ d^{3} \end{pmatrix}},$

where the objects differ in size. For instance, the first object may be un-sampled with dimension k and l and the second two objects may have a lower dimension than the first due to color sub-sampling. This, of course, only applies to color images. The quantization matrix can be seen as such a composition, too, i.e.,

$\overset{->}{Q} = {\begin{pmatrix} Q^{1} \\ Q^{2} \\ Q^{3} \end{pmatrix}.}$

Then, the quantization of the foregoing JPEG compression of the digital image data is inverted by multiplying the coefficient data d with the given quantization data Q, point-wise and for each component, and saving it again as coefficient data.

After this step the coefficient data

$\overset{->}{d} = \begin{pmatrix} d^{1} \\ d^{2} \\ d^{3} \end{pmatrix}$

is given as follows:

${d^{1} = \begin{pmatrix} {d_{1,1}^{1}Q_{1,1}^{1}} & \ldots & {d_{1,{{8\; l_{1}} - 1}}^{1}Q_{1,{{8\; l_{1}} - 1}}^{1}} \\ \vdots & \; & \vdots \\ {d_{{{8k_{1}} - 1},1}^{1}Q_{{{8k_{1}} - 1},1}^{1}} & \ldots & {d_{{{8k_{1}} - 1},{{8l_{1}} - 1}}^{1}Q_{{{8k_{1}} - 1},{{8l_{1}} - 1}}^{1}} \end{pmatrix}},{d^{2} = \begin{pmatrix} {d_{1,1}^{2}Q_{1,1}^{2}} & \ldots & {d_{1,{{8\; l_{2}} - 1}}^{2}Q_{1,{{8\; l_{2}} - 1}}^{2}} \\ \vdots & \; & \vdots \\ {d_{{{8k_{2}} - 1},1}^{2}Q_{{{8k_{2}} - 1},1}^{2}} & \ldots & {d_{{{8k_{2}} - 1},{{8l_{2}} - 1}}^{2}Q_{{{8k_{2}} - 1},{{8l_{2}} - 1}}^{2}} \end{pmatrix}},{d^{3} = {\begin{pmatrix} {d_{1,1}^{3}Q_{1,1}^{3}} & \ldots & {d_{1,{{8\; l_{3}} - 1}}^{3}Q_{1,{{8\; l_{3}} - 1}}^{3}} \\ \vdots & \; & \vdots \\ {d_{{{8k_{3}} - 1},1}^{3}Q_{{{8k_{3}} - 1},1}^{3}} & \ldots & {d_{{{8k_{3}} - 1},{{8l_{3}} - 1}}^{3}Q_{{{8k_{3}} - 1},{{8l_{3}} - 1}}^{3}} \end{pmatrix}.}}$

Pointwise, this step can be written as d_(i,j) ¹=d_(i,j) ¹Q_(i,j) ¹, d_(i,j) ²=d_(i,j) ²Q_(i,j) ² and d_(i,j) ³=d_(i,j) ³Q_(i,j) ³. The coefficient data d now contains the already dequantized data.

In a next step the image which will be used at the beginning of the following iterative process is initialized. In principle, any initialization of the image data {right arrow over (u)} is suitable, as for example {right arrow over (u)}=0. However, in order to already start with an image “close” to the desired output a standard-decompressed image is chosen. “Standard-decompression” here means image data that has been decompressed using conventional methods.

For that step, two operations are performed on the dequantized data {right arrow over (d)}, namely an inverse component-wise block-cosine transformation, denoted by IBDCT ({right arrow over (d)}) and an up-sampling operation, denoted by S⁻¹{right arrow over (d)} which is the pseudo-inversion of the color sub-sampling performed during the compression process. “Pseudo-inversion” here means that the inverse is not bijective, but only right-inverse.

To perform the inverse block-cosine transformation (IBDCT) on given coefficient data {right arrow over (d)} each image component is split into blocks of size 8×8 pixels and each of these blocks is then processed by an inverse cosine transformation. The resulting 8×8 blocks are then temporary stored in IBDCT({right arrow over (d)}), which is of the same data type as {right arrow over (d)}, at the position of the input block. For a given 8×8 image block (z_(i,j))_(0<i, j<8), the 8×8 block (IDCT(z)_(i,j))_(0<i, j<8) resulting from the IDCT operation is defined pointwise as

$\left( {{IDCT}(z)} \right)_{i,j} = {\sum\limits_{p,{q = 0}}^{7}{C_{p}C_{q}z_{p,q}{\cos \left( \frac{{\pi \left( {{2i} + 1} \right)}p}{16} \right)}{{\cos \left( \frac{{\pi \left( {{2j} + 1} \right)}q}{16} \right)}.}}}$

The constants C_(p) and C_(q) are defined as

$C_{p} = \left\{ \begin{matrix} \frac{1}{\sqrt{8}} & {{{{if}\mspace{14mu} p} = 0},} \\ \frac{1}{2} & {{{if}\mspace{14mu} 1} \leq p \leq 7} \end{matrix} \right.$

and analogously for C_(q).

This step is illustrated by an easy example, considering given data d₁ of size 16×16 written as

$d_{1} = {\begin{pmatrix} d_{0,0}^{1} & \ldots & d_{0,7}^{1} & d_{0,8}^{1} & \ldots & d_{0,15}^{1} \\ \vdots & \; & \vdots & \vdots & \; & \vdots \\ d_{7,0}^{1} & \ldots & d_{7,7}^{1} & d_{7,8}^{1} & \ldots & d_{7,15}^{1} \\ d_{8,0}^{1} & \ldots & d_{8,7}^{1} & d_{8,8}^{1} & \ldots & d_{8,15}^{1} \\ \vdots & \; & \vdots & \vdots & \; & \vdots \\ d_{15,0}^{1} & \ldots & d_{15,7}^{1} & d_{15,8}^{1} & \ldots & d_{15,15}^{1} \end{pmatrix}.}$

The inverse block-cosine transformation is then given as

${{IBDCT}\left( d_{1} \right)} = \begin{pmatrix} {{IDCT}\begin{pmatrix} d_{0,0}^{1} & \ldots & d_{0,7}^{1} \\ \vdots & \; & \vdots \\ d_{7,0}^{1} & \ldots & d_{7,7}^{1} \end{pmatrix}} & {{IDCT}\begin{pmatrix} d_{0,8}^{1} & \ldots & d_{0,15}^{1} \\ \vdots & \; & \vdots \\ d_{7,8}^{1} & \ldots & d_{7,15}^{1} \end{pmatrix}} \\ {{IDCT}\begin{pmatrix} d_{8,0}^{1} & \ldots & d_{8,7}^{1} \\ \vdots & \; & \vdots \\ d_{15,0}^{1} & \ldots & d_{15,7}^{1} \end{pmatrix}} & {{IDCT}\begin{pmatrix} d_{8,8}^{1} & \ldots & d_{8,15}^{1} \\ \vdots & \; & \vdots \\ d_{15,8}^{1} & \ldots & d_{15,15}^{1} \end{pmatrix}} \end{pmatrix}$

with the IDCT-operation performed as defined above.

For the up-sampling S⁻¹ it has to be considered that the components d¹, d² and d³ may have been sub-sampled and are saved in a lower resolution than the original image. Hence, the color components d¹, d² and d³ are potentially up-sampled, usually by repeating them up to the image-size.

Also, before the up-sampling process the color image was denoted as composition of three different image matrices while after the sub-sampling it is denoted as a vector-valued matrix. This results from the fact that before the up-sampling the three components may have different size and after up-sampling they all have the same size, making the vector-notation more convenient.

Generally, for given data

$\overset{\rightarrow}{d} = \begin{pmatrix} d^{1} \\ d^{2} \\ d^{3} \end{pmatrix}$

the pointwise assignment of the entries of the up-sampled data S⁻¹{right arrow over (d)} matrix is given by

$\left( {S^{- 1}d} \right)_{{{\frac{k}{k_{1}}i} + n},{{\frac{l}{l_{1}}j} + m}}^{1} = d_{ij}^{1}$ for i = 0  …  8k₁ − 1, j = 0  …  8l₁ − 1 and n = 0  …  (k/k₁) − 1, m = 0  …  (l/l₁) − 1; $\left( {S^{- 1}d} \right)_{{{\frac{k}{k_{2}}i} + n},{{\frac{l}{l_{2}}j} + m}}^{2} = d_{ij}^{2}$ for i = 0  …  8k₂ − 1, j = 0  …  8l₂ − 1 and n = 0  …  (k/k₂) − 1, m = 0  …  (l/l₂) − 1; $\left( {S^{- 1}d} \right)_{{{\frac{k}{k_{3}}i} + n},{{\frac{l}{l_{3}}j} + m}}^{3} = d_{ij}^{2}$ for i = 0  …  8k₃ − 1, j = 0  …  8l₃ − 1 and n = 0  …  (k/k₃) − 1, m = 0  …  (l/l₃) − 1.

In a forth step a number of variables is defined.

The extragradient image data {right arrow over (u)} is initialized by the image data {right arrow over (u)}, the scalars σ, τ are chosen arbitrarily such that the condition στ≦ 1/12 is satisfied and all other variables are initialized with 0 in all their components. It has to be noted that this is one exemplary initialization; in principle, any initialization is possible. For the extended image data {right arrow over (v)}, the extragradient extended image data {right arrow over (v)} and the two-vector valued matrix {right arrow over (p)} which are of two-vector matrix type, this can be written as

$\overset{\rightarrow}{v} = \left( {\begin{pmatrix} 0 \\ 0 \\ 0 \end{pmatrix}\begin{pmatrix} 0 \\ 0 \\ 0 \end{pmatrix}} \right)_{{0 \leq i < {8k}},{0 \leq j < {8l}}}$

and analogously for {right arrow over (v)} and {right arrow over (p)}.

The variable {right arrow over (q)} is of the three-vector matrix valued matrix type and looks like

$\overset{\rightarrow}{q} = {\left( {\begin{pmatrix} 0 \\ 0 \\ 0 \end{pmatrix}\begin{pmatrix} 0 \\ 0 \\ 0 \end{pmatrix}\begin{pmatrix} 0 \\ 0 \\ 0 \end{pmatrix}} \right)_{{0 \leq i < {8k}},{0 \leq j < {8l}}}.}$

The three-vector valued matrix q and the two-vector valued matrix p are the dual data that can be used as input for primal-dual methods.

The following operations A to G are repeated until a stop criterion (see below) is satisfied:

Step A: This step consists of a gradient evaluation 109, a summation/multiplication operation and a projection 110. It processes data stored in {right arrow over (p)}, {right arrow over (u)} and {right arrow over (v)}, the result is stored in {right arrow over (p)}.

The gradient evaluation V maps the extragradient image data {right arrow over (u)} to a two-vector valued matrix denoted by

${\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} = \left( {\begin{pmatrix} \left( {\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} \right)_{i,j}^{1,1} \\ \left( {\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} \right)_{i,j}^{2,1} \\ \left( {\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} \right)_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} \left( {\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} \right)_{i,j}^{1,2} \\ \left( {\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} \right)_{i,j}^{2,2} \\ \left( {\bigtriangledown \overset{\rightarrow}{\overset{\_}{u}}} \right)_{i,j}^{3,2} \end{pmatrix}} \right)_{{0 \leq i < {8k}},{0 \leq j < {8l}}}$

where the entries are defined as follows:

( V {right arrow over (u)})_(i,j) ^(n,1)=(δ_(x)+ u ^(n))_(i,j),

( V {right arrow over (u)}) _(i,j) ^(n,2)=(δ_(y+) ū ^(n))_(i,j)

for n=1, 2, 3. The operators δ_(x+) and δ_(y+), which are defined on a scalar valued matrix z=(z_(i,j))_(0≦i<8 k, 0≦j<8 l) result again in a scalar valued matrix ((δ_(x+)z)_(i,j))_(0≦i<8 k, 0≦j<8 l) and ((δ_(y+)z)_(i,j))_(0≦i<8 k, 0≦j<8 l), respectively, which is given pointwise as

$\left( {\delta_{x +}z} \right)_{i,j} = \left\{ {{\begin{matrix} \left( {z_{{i + 1},j} - z_{i,j}} \right) & {{{{if}\mspace{14mu} 0} \leq i < {{8k} - 1}},} \\ 0 & {{{{if}\mspace{14mu} i} = {{8k} - 1}},} \end{matrix}{{and}\left( {\delta_{y +}z} \right)}_{i,j}} = \left\{ \begin{matrix} \left( {z_{i,{j + 1}} - z_{i,j}} \right) & {{{{if}\mspace{14mu} 0} \leq j < {{8l} - 1}},} \\ 0 & {{{if}\mspace{14mu} j} = {{8l} - 1.}} \end{matrix} \right.} \right.$

The summation/multiplication operation uses this data and yields a two-vector valued matrix which is stored in p, with its entries being defined pointwise by

p _(i,j) ^(m,n) =p _(i,j) ^(m,n)+σ((∇ {right arrow over (u)} )_(i,j) ^(m,n) − v _(i,j) ^(m,n)).

The projection step then modifies each two-vector entry of p, denoted by

${{\overset{\rightarrow}{p}}_{i,j} = \left( {\begin{pmatrix} (p)_{i,j}^{1,1} \\ (p)_{i,j}^{2,1} \\ (p)_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} (p)_{i,j}^{1,2} \\ (p)_{i,j}^{2,2} \\ (p)_{i,j}^{3,2} \end{pmatrix}} \right)},$

in a way that its norm, denoted by

${{n\; 1_{i,j}} = {\left( {\begin{pmatrix} (p)_{i,j}^{1,1} \\ (p)_{i,j}^{2,1} \\ (p)_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} (p)_{i,j}^{1,2} \\ (p)_{i,j}^{2,2} \\ (p)_{i,j}^{3,2} \end{pmatrix}} \right)}},$

is less or equal a given α₁. As a result, the new entries for the two-vector valued matrix p are given by

$p_{i,j}^{m,n} = {\frac{p_{i,j}^{m,n}}{\max \left\{ {1,\frac{n\; 1_{i,j}}{\alpha_{1}}} \right\}}.}$

A possible choice for the norm n1_(i,j) of the entry {right arrow over (p)}_(i,j) is

${n\;}_{i,j} = {\sqrt{\left( p_{i,j}^{1,1} \right)^{2} + \left( p_{i,j}^{2,1} \right)^{2} + \left( p_{i,j}^{3,1} \right)^{2} + \left( p_{i,j}^{1,2} \right)^{2} + \left( p_{i,j}^{2,2} \right)^{2} + \left( p_{i,j}^{3,2} \right)^{2}}.}$

Step B: This step is similar to step A. It consists of a higher dimensional gradient evaluation 109, a summation/multiplication operation and a projection operation 110.

The gradient evaluations maps the ex tragradient extended image data {right arrow over (v)} in the form of a two-vector valued matrix to a three-vector valued matrix ε( {right arrow over (v)}) which is defined pointwise by

$\left( {ɛ\overset{\rightarrow}{\overset{\_}{v}}} \right)_{i,j} = {\left( {\begin{pmatrix} \left( {\delta_{x -}{\overset{\_}{v}}^{1,1}} \right)_{i,j} \\ \left( {\delta_{x -}{\overset{\_}{v}}^{2,1}} \right)_{i,j} \\ \left( {\delta_{x -}{\overset{\_}{v}}^{3,1}} \right)_{i,j} \end{pmatrix},\begin{pmatrix} \left( {\delta_{y -}{\overset{\_}{v}}^{1,2}} \right)_{i,j} \\ \left( {\delta_{y -}{\overset{\_}{v}}^{2,2}} \right)_{i,j} \\ \left( {\delta_{y -}{\overset{\_}{v}}^{3,2}} \right)_{i,j} \end{pmatrix},\begin{pmatrix} \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{1,1}} + {\delta_{x -}{\overset{\_}{v}}^{1,2}}}{2} \right)_{i,j} \\ \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{2,1}} + {\delta_{x -}{\overset{\_}{v}}^{2,2}}}{2} \right)_{i,j} \\ \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{3,1}} + {\delta_{x -}{\overset{\_}{v}}^{3,2}}}{2} \right)_{i,j} \end{pmatrix}} \right).}$

The operators δ_(x−), δ_(y−) are also defined on scalar valued matrices z=(z_(i,j))_(0≦i<8 k, 0≦j<8 l) as

$\left( {\delta_{x -}z} \right)_{i,j} = \left\{ {{\begin{matrix} {- z_{{i - 1},j}} & {{{if}\mspace{14mu} i} = {{8k} - 1}} \\ \left( {z_{i,j} - z_{{i - 1},j}} \right) & {{{{if}\mspace{14mu} 0} < i < {{8k} - 1}},} \\ z_{i,j} & {{{if}\mspace{14mu} i} = 0} \end{matrix}\left( {\delta_{y -}z} \right)_{i,j}} = \left\{ \begin{matrix} {- z_{i,{j - 1}}} & {{{if}\mspace{14mu} i} = {{8l} - 1}} \\ \left( {z_{i,j} - z_{i,{j - 1}}} \right) & {{{if}\mspace{14mu} 0} < j < {{8l} - 1}} \\ z_{i,j} & {{{if}\mspace{14mu} j} = 0.} \end{matrix} \right.} \right.$

The summation/multiplication operation uses this data and yields a three-vector valued matrix which is stored in {right arrow over (q)} with its entries being defined pointwise by

q _(i,j) ^(m,n) =q _(i,j) ^(m,n)+σ(ε {right arrow over (v)})_(i,j) ^(m,n).

The projection step then modifies each three-vector entry of the three-vector valued matrix {right arrow over (q)}, denoted by

${\overset{\rightarrow}{q}}_{i,j} = \left( {\begin{pmatrix} (q)_{i,j}^{1,1} \\ (q)_{i,j}^{2,1} \\ (q)_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} (q)_{i,j}^{1,2} \\ (q)_{i,j}^{2,2} \\ (q)_{i,j}^{3,2} \end{pmatrix},\begin{pmatrix} (q)_{i,j}^{1,3} \\ (q)_{i,j}^{2,3} \\ (q)_{i,j}^{3,3} \end{pmatrix}} \right)$

in a way that its norm, denoted by

${n\; 2_{i,j}} = {\left( {\begin{pmatrix} (q)_{i,j}^{1,1} \\ (q)_{i,j}^{2,1} \\ (q)_{i,j}^{3,1} \end{pmatrix},\begin{pmatrix} (q)_{i,j}^{1,2} \\ (q)_{i,j}^{2,2} \\ (q)_{i,j}^{3,2} \end{pmatrix},\begin{pmatrix} (q)_{i,j}^{1,3} \\ (q)_{i,j}^{2,3} \\ (q)_{i,j}^{3,3} \end{pmatrix}} \right)}$

is less or equal a given α₀. As result, the new entry is given by

$q_{i,j}^{m,n} = {\frac{q_{i,j}^{m,n}}{\max \left\{ {1,\frac{n\; 2_{i,j}}{\alpha_{0}}} \right\}}.}$

A possible choice of the norm n2_(i,j) of the entry {right arrow over (q)}_(i,j) is

${n\; 2_{i,j}} = {\sqrt{\begin{matrix} {\left( q_{i,j}^{1,1} \right)^{2} + \left( q_{i,j}^{2,1} \right)^{2} + \left( q_{i,j}^{3,1} \right)^{2} + \left( q_{i,j}^{1,2} \right)^{2} +} \\ {\left( q_{i,j}^{2,2} \right)^{2} + \left( q_{i,j}^{3,2} \right)^{2} + {2\left( {\left( q_{i,j}^{1,3} \right)^{2} + \left( q_{i,j}^{2,3} \right)^{2} + \left( q_{i,j}^{3,3} \right)^{2}} \right)}} \end{matrix}}.}$

Step C: This step evaluates a divergence operator 111, performs a summation/multiplication operation and stores the result in a first new set of image data {right arrow over (u)}_(new). The divergence operator uses the two-vector valued matrix {right arrow over (p)} as input and results in a vector valued matrix div {right arrow over (p)} which is defined pointwise by

(div {right arrow over (p)})_(i,j) ^(n)=(δ_(x−) p ^(n,1))_(i,j)+(δ_(y−) p ^(n,2))_(i,j),

where δ_(x−), δ_(y−) is defined as above. The summation/multiplication step then results in a vector-valued matrix which is stored in {right arrow over (u)}_(new), defined pointwise by

(u _(new))_(i,j) ^(n) =u _(i,j) ^(n)+τ(div {right arrow over (p)})_(i,j) ^(n).

Step D: This step is similar to step C. It evaluates a divergence operator 111, performs a summation/multiplication operation and stores the result in a second new set of image data {right arrow over (v)}_(new). The divergence operator uses the three-vector valued matrix {right arrow over (q)} as input and results in a two-vector-valued matrix div {right arrow over (q)}, which is defined pointwise by

(div {right arrow over (q)})_(i,j) ^(n,1)=(δ_(x−) q ^(n,1))_(i,j)+(δ_(y+) q ^(n,3))_(i,j),

(div {right arrow over (q)})_(i,j) ^(n,2)=(δ_(x+) q ^(n,3))_(i,j)′(δ_(y+) q ^(n,2))_(i,j)

where δ_(x+), δ_(y+) are defined as above. The summation/multiplication step then results in a two vector-valued matrix, which is stored in {right arrow over (v)}_(new), defined pointwise by

(v _(new))_(i,j) ^(m,n) =v _(i,j) ^(m,n)+τ(p+div {right arrow over (q)})_(i,j) ^(m,n).

Step E: This step describes a project-to-data operation which uses the vector-valued matrix {right arrow over (u)}_(new) and writes the result proj_(D)({right arrow over (u)}_(new)) again to {right arrow over (u)}_(new). In order to get the result of this projection-operation, several computations have to be performed. At first, a sub-sampling operation 112 is performed on the first new set of image data {right arrow over (u)}_(new) which results in three, temporary stored, scalar-valued matrices (S{right arrow over (u)}_(new))¹, (S{right arrow over (u)}_(new))², (S{right arrow over (u)}_(new))³. The entries of these matrices depend on the sub-sampled image size parameters k₁, k₂, k₃, l₁, l₂, l₃ and are defined by

${\left( {S{\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{1} = {\frac{k_{1}l_{1}}{kl}{\sum\limits_{m = {i\frac{k}{k_{1}}}}^{{{({i + 1})}\frac{k}{k_{1}}} - 1}\; {\sum\limits_{n = {j\frac{l}{l_{1}}}}^{{{({j + 1})}\frac{l}{l_{1}}} - 1}\; \left( u_{new} \right)_{n,m}^{1}}}}},{\left( {S{\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{2} = {\frac{k_{2}l_{2}}{kl}{\sum\limits_{m = {i\frac{k}{k_{2}}}}^{{{({i + 1})}\frac{k}{k_{2}}} - 1}\; {\sum\limits_{n = {i\frac{l}{l_{2}}}}^{{{({i + 1})}\frac{l}{l_{2}}} - 1}\; \left( u_{new} \right)_{n,m}^{2}}}}},{\left( {S{\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{3} = {\frac{k_{3}l_{3}}{kl}{\sum\limits_{m = {i\frac{k}{k_{3}}}}^{{{({i + 1})}\frac{k}{k_{3}}} - 1}\; {\sum\limits_{n = {i\frac{l}{l_{3}}}}^{{{({i + 1})}\frac{l}{l_{3}}} - 1}\; {\left( u_{new} \right)_{n,m}^{3}.}}}}}$

Using these matrices as input data, a component-wise block-wise discrete cosine transformation 113 is performed and temporary stored in three scalar-valued matrices (BDCT(S{right arrow over (u)}_(new)))¹, (BDCT(S{right arrow over (u)}_(new)))² and (BDCT(S{right arrow over (u)}_(new))³. This is done by first splitting each scalar-valued input matrix into blocks of size 8×8-entries, then applying a discrete cosine transformation (DCT) on each of these 8×8 blocks and storing the resulting 8×8 block at the same position as the input block in (BDCT(S{right arrow over (u)}_(new)))^(n). For a 8×8-block (z_(i,j))_(0≦i, j<7) the 8×8-block DCT(z_(i,j)))_(0≦i, j<7) resulting from the discrete cosine transformation is defined pointwise as

${{{DCT}(z)}_{i,j} = {C_{i}C_{j}{\sum\limits_{n,{m = 0}}^{7}\; {z_{n,m}{\cos \left( \frac{{\pi \left( {{2\; n} + 1} \right)}i}{16} \right)}{\cos \left( \frac{{\pi \left( {{2\; m} + 1} \right)}j}{16} \right)}}}}},{where}$ $C_{s} = \left\{ {{{\begin{matrix} \frac{1}{\sqrt{8}} & {{{{if}\mspace{14mu} s} = 0},} \\ \frac{1}{2} & {{{if}\mspace{14mu} 1} \leq s \leq 7} \end{matrix}{for}\mspace{14mu} s} = i},{j.}} \right.$

The scalar-valued matrices (BDCT(S{right arrow over (u)}_(new)))^(n), n=1, 2, 3, which are denoted as

${{\overset{\rightarrow}{d}\; t} = \begin{pmatrix} {dt}^{1} \\ {dt}^{2} \\ {dt}^{3} \end{pmatrix}},$

are then processed in a projection operation 114 which results again in three scalar-valued matrices

${{proj}_{d,Q}\left( {\overset{\rightarrow}{d}\; t} \right)}_{i,j}^{n} = {{dt}_{i,j}^{n} - {\max \left\{ {0,{{dt}_{i,j}^{n} - d_{i,j}^{n} - \frac{Q_{i,j}^{n}}{2}}} \right\}} + {\max {\left\{ {0,{d_{i,j}^{n} - \frac{Q_{i,j}^{n}}{2} - {dt}_{i,j}^{n}}} \right\}.}}}$

Using results of the operations IBDCT 115 and 116 as previously defined in a reassembling operation 117, the result of the projection-to-data operation , denoted by proj_(D)({right arrow over (u)}_(new)), is then given pointwise as

proj_(D)({right arrow over (u)} _(new))_(i,j) ^(n)=(u _(new))_(i,j) ^(n)+(S ⁻¹(IBDCT(proj_(d,Q)({right arrow over (d)}t))−S{right arrow over (u)} _(new) ))_(i,j) ^(n) and stored in u _(new).

Step F: The next operation is an extragradient-update step 118. It processes data stored in {right arrow over (u)}, {right arrow over (v)}, {right arrow over (u)}_(new), {right arrow over (v)}_(new) and saves the result in {right arrow over (u)}, {right arrow over (v)} which are then given pointwise as

ū _(i,j) ^(n)=2(u _(new))_(i,j) ^(n) −u _(i,j) ^(n),

v _(i,j) ^(m,n)=2(v _(new))_(i,j) ^(m,n) −v _(i,j) ^(m,n).

Step G: In this step, the “old” image data it and extended image data {right arrow over (v)} is overwritten by the first new set of image data {right arrow over (u)}_(new) and the second new set of image data {right arrow over (v)}_(new) which is done pointwise by

u _(i,j) ^(n)=(u _(new))_(i,j) ^(n),

v _(i,j) ^(m,n)=(v _(new))_(i,j) ^(m,n).

Step H: After each iteration of steps A-G a previously defined stopping criterion is tested 119. If the stopping criterion is fulfilled the iterations are halted; if not the iterations are continued and the process starts anew at A.

Possible choices for the stopping criterion are a certain maximal number of iterations, a certain maximal amount of time consumed by the process or a data-dependent criterion like an ε-threshold or the like.

If the stopping criterion is fulfilled the process 120 returns the image saved in {right arrow over (u)}_(new) as de-compression of the initial digital image data, i.e., the given JPEG-object. With that, the process is stopped. 

1. Method for the reconstruction of compressed digital image data which has been lossy-compressed by transform coding and quantization, comprising the steps of a) Inputting JPEG-compressed digital image data (J_(comp)); b) Iteratively performing updates of an extended image model, thereby reducing the total generalized variation functional constrained to all images matching the transform coded, quantized data; c) Outputting a decompressed digital image (u).
 2. The method of claim 1, wherein in step b) the extended image model comprises an image, an extended image, a first dual variable and a second dual variable, an extragradient image and an extragradient extended image.
 3. The method of claim 1, wherein step b) comprises the following sub-steps: b0) Initializing image data (u), extended image data (v), extragradient image data (ū), extragradient extended image data ( v), a two-vector valued matrix (p) and a three-vector valued matrix (q) which are used as a starting point for the iteration procedure in steps b1) to b8); b1) Applying a gradient evaluation and a summation/multiplication operation on the extragradient image data (ū), the extragradient extended image data ( v) and the two-vector valued matrix (p), thereby updating the two-vector valued matrix (p) and then applying a projection operation on said matrix (p) in a way that its norm (n1_(i,j)) is less or equal to a first parameter (α₁) and each two-vector entry of the two-vector valued matrix (p) is modified accordingly; b2) Applying a higher dimensional gradient evaluation and a summation/multiplication operation on the extragradient extended image data ( v) and the three-vector valued matrix (q), thereby updating the results in the three-vector valued matrix (q) and then applying a projection operation on said matrix (q) in a way that its norm (n2_(i,j)) is less or equal to a second parameter (α₀) and each three-vector entry of the three-vector valued matrix (q) is modified accordingly; b3) Applying a divergence operator and a summation/multiplication operation on the image data (u) and the two-vector valued matrix (p) of step b1), storing the data in a first new set of image data (u_(new)); b4) Applying a divergence operator and a summation/multiplication operation on the extended image data (v), the two-vector valued matrix (p) of step b1) and the three-vector-valued matrix (q) of step b2), storing the data in a second new set of image data (v_(new)); b5) Updating the first new set of image data (u_(new)) pf step b3) by applying a sub-sampling operation S, a component-wise block-wise discrete cosine transformation and a projection operation on the first new set of image data, followed by an inverse component-wise block-cosine transformation IBDCT, an up-sampling operation S⁻¹ and a summation operation; b6) Updating the extragradient image data (ū) and the extragradient extended image data ( v) by an extragradient-update step using the image data (u), extended image data (v), first (u_(new)) and second new set of image data (v_(new)); b7) Overwriting image data (u) and extended image data (v) with of the first (u_(new)) and second new set of image data (v_(new)), respectively, resulting from the steps b1) to b6); b8) Evaluating a stopping-criterion; b9) Reiterating steps b1) to b8) if stopping criterion is not fulfilled or proceeding to step c) if stopping criterion is fulfilled.
 4. The method of claim 3, wherein in step b0) a standard-decompressed image is selected, performing an inverse component-wise block-cosine transformation IBDCT and an up-sampling operation S⁻¹ on the dequantized data (d).
 5. Method of claim 3, wherein in step b1) the gradient-evaluation maps extragradient image data ( {right arrow over (u)}) to a two-vector-valued matrix with the entries being defined as follows:) ( V {right arrow over (u)} )_(i,j) ^(n,1)=(δ_(x+) ū ^(n))_(i,j), and ( V {right arrow over (u)} )_(i,j) ^(n,2)=(δ_(y+) ū ^(n))_(i,j), with 0≦i<8 k and 0≦j<8l, k and l corresponding to the vertical and horizontal number of 8×8 data blocks and n=1, 2, 3, the operators each resulting in a scalar-valued matrix with $\left( {\delta_{x +}\overset{\_}{u}} \right)_{i,j} = \left\{ {{\begin{matrix} \left( {{\overset{\_}{u}}_{{i + 1},j} - {\overset{\_}{u}}_{i,j}} \right) & {{{{if}\mspace{14mu} 0} \leq i < {{8\; k} - 1}},} \\ 0 & {{{{if}\mspace{14mu} i} = {{8\; k} - 1}},} \end{matrix}{{and}\left( {\delta_{y +}\overset{\_}{u}} \right)}_{i,j}} = \left\{ {\begin{matrix} \left( {{\overset{\_}{u}}_{i,{j + 1}} - {\overset{\_}{u}}_{i,j}} \right) & {{{{if}\mspace{14mu} 0} \leq j < {{8\; l} - 1}},} \\ 0 & {{{if}\mspace{14mu} j} = {{8\; l} - 1}} \end{matrix}{and}} \right.} \right.$ the summation/multiplication operation yields a two-vector valued matrix p_(i,j) ^(m,n)) with p_(i,j) ^(m,n)=p_(i,j) ^(m,n)+σ(( V {right arrow over (u)})_(i,j) ^(m,n)− v _(i,j) ^(m,n)), σ being a first scalar and v _(i,j) ^(m,n) being the extragradient extended image data.
 6. Method of claim 3, wherein in step b1) the modification of the two-vector valued matrix (p_(i,j) ^(m,n)) is calculated with ${p_{i,j}^{m,n} = \frac{p_{i,j}^{m,n}}{\max \left\{ {1,\frac{n\; 1_{i,j}}{\alpha_{1}}} \right\}}},$ α₁ being a first parameter and n1 being a norm.
 7. Method of claim 6, wherein the norm (n1_(i,j)) is calculated as n1_(i,j)=√{square root over ((p_(i,j) ^(1,1))²+(p_(i,j) ^(2,1))²+(p_(i,j) ^(3,1))²+(p_(i,j) ^(1,2))²+(p_(i,j) ^(2,2))²+(p_(i,j) ^(3,2))²)}{square root over ((p_(i,j) ^(1,1))²+(p_(i,j) ^(2,1))²+(p_(i,j) ^(3,1))²+(p_(i,j) ^(1,2))²+(p_(i,j) ^(2,2))²+(p_(i,j) ^(3,2))²)}{square root over ((p_(i,j) ^(1,1))²+(p_(i,j) ^(2,1))²+(p_(i,j) ^(3,1))²+(p_(i,j) ^(1,2))²+(p_(i,j) ^(2,2))²+(p_(i,j) ^(3,2))²)}{square root over ((p_(i,j) ^(1,1))²+(p_(i,j) ^(2,1))²+(p_(i,j) ^(3,1))²+(p_(i,j) ^(1,2))²+(p_(i,j) ^(2,2))²+(p_(i,j) ^(3,2))²)}{square root over ((p_(i,j) ^(1,1))²+(p_(i,j) ^(2,1))²+(p_(i,j) ^(3,1))²+(p_(i,j) ^(1,2))²+(p_(i,j) ^(2,2))²+(p_(i,j) ^(3,2))²)}{square root over ((p_(i,j) ^(1,1))²+(p_(i,j) ^(2,1))²+(p_(i,j) ^(3,1))²+(p_(i,j) ^(1,2))²+(p_(i,j) ^(2,2))²+(p_(i,j) ^(3,2))²)}. ${n\; 1_{i,j}} = {\sqrt{\left( p_{i,j}^{1,1} \right)^{2} + \left( p_{i,j}^{2,1} \right)^{2} + \left( p_{i,j}^{3,1} \right)^{2} + \left( p_{i,j}^{1,2} \right)^{2} + \left( p_{i,j}^{2,2} \right)^{2} + \left( p_{i,j}^{3,2} \right)^{2}}.}$
 8. Method of claim 3, wherein in step b2) the higher dimensional gradient evaluation maps extragradient extended image data ( v) to a three-vector valued matrix being defined as follows: ${\left( {ɛ\; \overset{\rightarrow}{\overset{\_}{v}}} \right)_{i,j} = \left( {\begin{pmatrix} \left( {\delta_{x -}{\overset{\_}{v}}^{1,1}} \right)_{i,j} \\ \left( {\delta_{x -}{\overset{\_}{v}}^{2,1}} \right)_{i,j} \\ \left( {\delta_{x -}{\overset{\_}{v}}^{3,1}} \right)_{i,j} \end{pmatrix},\begin{pmatrix} \left( {\delta_{y -}{\overset{\_}{v}}^{1,2}} \right)_{i,j} \\ \left( {\delta_{y -}{\overset{\_}{v}}^{2,2}} \right)_{i,j} \\ \left( {\delta_{y -}{\overset{\_}{v}}^{3,2}} \right)_{i,j} \end{pmatrix},\begin{pmatrix} \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{1,1}} + {\delta_{x -}{\overset{\_}{v}}^{1,2}}}{2} \right)_{i,j} \\ \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{2,1}} + {\delta_{x -}{\overset{\_}{v}}^{2,2}}}{2} \right)_{i,j} \\ \left( \frac{{\delta_{y -}{\overset{\_}{v}}^{3,1}} + {\delta_{x -}{\overset{\_}{v}}^{3,2}}}{2} \right)_{i,j} \end{pmatrix}} \right)},$ the operators being defined as $\left( {\delta_{x -}\overset{\_}{v}} \right)_{i,j} = \left\{ {{\begin{matrix} {- {\overset{\_}{v}}_{{i - 1},j}} & {{{if}\mspace{14mu} i} = {{8\; k} - 1}} \\ \left( {{\overset{\_}{v}}_{i,j} - {\overset{\_}{v}}_{{i - 1},j}} \right) & {{{if}\mspace{14mu} 0} < i < {{8\; k} - 1}} \\ {\overset{\_}{v}}_{i,j} & {{{if}\mspace{14mu} i} = 0} \end{matrix}{{and}\left( {\delta_{y -}\overset{\_}{v}} \right)}_{i,j}} = \left\{ {\begin{matrix} {- {\overset{\_}{v}}_{i,{j - 1}}} & {{{if}\mspace{14mu} i} = {{8\; l} - 1}} \\ \left( {{\overset{\_}{v}}_{i,j} - {\overset{\_}{v}}_{i,{j - 1}}} \right) & {{{if}\mspace{14mu} 0} < j < {{8\; l} - 1}} \\ {\overset{\_}{v}}_{i,j} & {{{if}\mspace{14mu} j} = 0.} \end{matrix},} \right.} \right.$ and the summation/multiplication operation yields a three-vector valued matrix (q_(i,j) ^(m,n)) with q_(i,j) ^(m,n)=q_(i,j) ^(m,n)+σ(ε {right arrow over (v)})_(i,j) ^(m,n).
 9. Method of claim 3, wherein in step b2) the modification of the three-vector valued matrix (q_(i,j) ^(m,n)) is calculated with ${q_{i,j}^{m,n} = \frac{q_{i,j}^{m,n}}{\max \left\{ {1,\frac{n\; 2_{i,j}}{\alpha_{0}}} \right\}}},$ α₀ being a second parameter and n2 being a norm.
 10. Method of claim 9, wherein the norm (n2_(i,j)) is calculated with ${n\; 2_{i,j}} = {\sqrt{\begin{matrix} {\left( q_{i,j}^{1,1} \right)^{2} + \left( q_{i,j}^{2,1} \right)^{2} + \left( q_{i,j}^{3,1} \right)^{2} + \left( q_{i,j}^{1,2} \right)^{2} +} \\ {\left( q_{i,j}^{2,2} \right)^{2} + \left( q_{i,j}^{3,2} \right)^{2} + {2\left( {\left( q_{i,j}^{1,3} \right)^{2} + \left( q_{i,j}^{2,3} \right)^{2} + \left( q_{i,j}^{3,3} \right)^{2}} \right)}} \end{matrix}}.}$
 11. Method of claim 3, wherein in step b3) a divergence operator is applied to the two-vector valued matrix (p_(i,j) ^(m,n)) of step b1) as follows: (div {right arrow over (p)})_(i,j) ^(n)=(δ_(x−)p^(n,1))_(i,j)+(δ_(y−)p^(n,2))_(i,j), and the summation/multiplication operation yields a first new set of image data (ū_(new)) as a vector-valued matrix with (u_(new))_(i,j) ^(n)=u_(i,j) ^(n)+τ(div {right arrow over (p)})_(i,j) ^(n), τ being a second scalar.
 12. Method of claim 3, wherein in step b4) a divergence operator is applied to the three-vector valued matrix (q_(i,j) ^(m,n)) of step b2) as follows: (div {right arrow over (q)})_(i,j) ^(n,1)=(δ_(x+) q ^(n,1))_(i,j)+(δ_(y+) q ^(n,3))_(i,j), (div {right arrow over (q)})_(i,j) ^(n,2)=(δ_(x+)q^(n,3))_(i,j)+(δ_(y+)q^(n,2))_(i,j) and the summation/multiplication operation yields a second new set of image data ({right arrow over (v)}_(new)) as a two-vector valued matrix with (v_(new))_(i,j) ^(m,n)=v_(i,j) ^(m,n)+τ({right arrow over (p)}+div {right arrow over (q)})_(i,j) ^(m,n).
 13. Method of claim 5, wherein the first (σ) and the second scalar (τ) satisfy the condition στ≦ 1/12.
 14. Method of claim 3, wherein in step b5) a sub-sampling operation is applied to the first new set of image data ({right arrow over (u)}_(new)) resulting in scalar-valued matrices as follows: ${\left( {S{\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{1} = {\frac{k_{1}l_{1}}{kl}{\sum\limits_{m = {i\frac{k}{k_{1}}}}^{{{({i + 1})}\frac{k}{k_{1}}} - 1}\; {\sum\limits_{n = {j\frac{l}{l_{1}}}}^{{{({j + 1})}\frac{l}{l_{1}}} - 1}\; \left( u_{new} \right)_{n,m}^{1}}}}},,$ $\left( {S{\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{2} = {\frac{k_{2}l_{2}}{kl}{\sum\limits_{m = {i\frac{k}{k_{2}}}}^{{{({i + 1})}\frac{k}{k_{2}}} - 1}\; {\sum\limits_{n = {i\frac{l}{l_{2}}}}^{{{({i + 1})}\frac{l}{l_{2}}} - 1}\; \left( u_{new} \right)_{n,m}^{2}}}}$ ${\left( {S{\overset{\rightarrow}{u}}_{new}} \right)_{i,j}^{3} = {\frac{k_{3}l_{3}}{kl}{\sum\limits_{m = {i\frac{k}{k_{3}}}}^{{{({i + 1})}\frac{k}{k_{3}}} - 1}\; {\sum\limits_{n = {i\frac{l}{l_{3}}}}^{{{({i + 1})}\frac{l}{l_{3}}} - 1}\; \left( u_{new} \right)_{n,m}^{3}}}}},$ with k₁, k₂, k₃, l₁, l₂, l₃ being the sub-sampled image size parameters, then a component-wise block-wise discrete cosine transformation is performed resulting in scalar-valued matrices (BDCT(S{right arrow over (u)}_(new)))¹, (BDCT(S{right arrow over (u)}_(new)))² (BDCT(S{right arrow over (u)}_(new)))³ and which are then processed in a projection operation with proj_(d,Q)(BDCT(S{right arrow over (u)} _(new)))_(i,j) ^(n)=BDCT(S{right arrow over (u)} _(new))_(i,j) ^(n)−max {0, BDCT(S{right arrow over (u)} _(new))_(i,j) ^(n)−(d _(i,j) ^(n)½)Q _(i,j) ^(n)}++max {0, d _(i,j) ^(n)−½)Q _(i,j) ^(n)−BDCT(S{right arrow over (u)} _(new))_(i,j) ^(n)} for n=1, 2, 3, resulting in three scalar-valued matrices which are then inversely transformed IBDCT and inversely sub-sampled S⁻¹, resulting in proj_(D)({right arrow over (u)} _(new))_(i,j) ^(n)=(u _(new))_(i,j) ^(n)+(S ⁻¹(IBDCT(proj_(d,Q)((S{right arrow over (u)} _(new))))−S{right arrow over (u)} _(new)))_(i,j) ^(n) which is then stored as first new set of image data ({right arrow over (u)}_(new)).
 15. Method of claim 3, wherein the extragradient update-step is performed by processing image data ({right arrow over (u)}, {right arrow over (v)}) and the first ({right arrow over (u)}_(new)) and second new set of image data ({right arrow over (v)}_(new)) as follows: ū _(i,j) ^(n)=2(u _(new))_(i,j) ^(n) −u _(i,j) ^(n), v _(i,j) ^(m,n)=2(v _(new))_(i,j) ^(m,n) −v _(i,j) ^(m,n) resulting in new sets of extragradient image data ( {right arrow over (u)}) and extragradient extended image data ( {right arrow over (v)}).
 16. Method of claim 1, wherein the JPEG-compressed digital image data (J_(comp)) is a greyscale image consisting of a matrix (u_(i,j)), each matrix element indicating the intensity of the respective pixel, or a color image in an RGB-color-space-setting consisting of a three-vector valued matrix $\left( {\left( {\overset{\rightarrow}{u}}_{i,j} \right) = \begin{pmatrix} u_{i,j}^{1} \\ u_{i,j}^{2} \\ u_{i,j}^{3} \end{pmatrix}} \right),$ each matrix vector indicating the amount of Red, Green and Blue light in the corresponding (i, j)-th pixel of the image, or a color image in an YCbCr-color-space-setting consisting of a three-vector valued matrix $\left( {\left( {\overset{\rightarrow}{u}}_{i,j} \right) = \begin{pmatrix} u_{i,j}^{1} \\ u_{i,j}^{2} \\ u_{i,j}^{3} \end{pmatrix}} \right),$ the matrix vectors indicating the brightness, the blue-difference and the red-difference of the corresponding pixel. 